#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <signal.h>
#include "mpi.h"

double f(double a){
    return(4*sqrt(1-a*a));
}
//extern f(double a);
double mid_point_rule(double start, double delta_x){

    return delta_x*f(start+(delta_x/2.))  ;
}
double rectangular_rule(double start, double delta_x){

    return delta_x*f(start)    ;
}
double simpsons_rule(double start, double delta_x){

    return(delta_x/6.0)*(f(start)+4*f(start+(delta_x/2.))+f(start+delta_x)) ;
}
double (*rule)(double start, double delta_x);

int main(int argc,char *argv[],char *envp[]){
    int done = 0, n, myid, numprocs, mpi_error,namelen,i_counter;
    int buffer[2],local_intervals;
    double local_start;
    double PI_to_25dp = 3.141592653589793238462643;
    double mypi, pi, delta_x;
    double startwtime, endwtime,beginwtime;
    char processor_name[MPI_MAX_PROCESSOR_NAME];

    MPI_Init(&argc,&argv);
    MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
    MPI_Comm_rank(MPI_COMM_WORLD,&myid);
    MPI_Get_processor_name(processor_name,&namelen);

    printf("Process number %d on processor %s \n", myid, processor_name, numprocs);
    fflush(stdout);
    n = 0;
    if (myid == 0) {
        buffer[0]=atoi(getenv("mpiMethod"));
        buffer[1]=atoi(getenv("mpiCount"));
        startwtime = MPI_Wtime();        
    }//end setup by the master process..0
    MPI_Bcast(buffer, 2, MPI_INT, 0, MPI_COMM_WORLD);
    if (buffer[0] == 0)
        done = 1;
    else {
        n=buffer[1];
        switch (buffer[0]) {
        case 0:
            done=1;
            break;
        case 1:
            rule=&rectangular_rule;
            break;
        case 2:
            rule=&mid_point_rule;
            break;
        case 3:
            rule=&simpsons_rule;
            break;
        }
        /*******set up this processes part of the calculation********/
        delta_x   = 1.0 / (double) n;
        local_start=(1.0/numprocs)*(myid);
        local_intervals=n/numprocs;
        /******* The next statments are just for diagnostics ********/
        //printf("process=%d local start = %.16f local intervals = %d\n",myid,local_start,local_intervals);
        //fflush(stdout);  
        /**********************************************************/
        mypi=0.0;
        i_counter=0;
        for (i_counter=0;i_counter<local_intervals;i_counter++)
            mypi+=( *rule)(local_start+i_counter*delta_x,delta_x);
        /*  The next statment---either sends mypi back to 0 or if 0 does the calculation*/
        MPI_Reduce(&mypi, &pi, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);

        if (myid == 0) {
            printf("\npi is approximately %.16f, Relative Error is %16.8e\n",
                   pi, (double)100 * (pi - PI_to_25dp)/PI_to_25dp);
            fflush(stdout);
            endwtime = MPI_Wtime();
            printf("wall clock time = %f\n", endwtime-startwtime);
            fflush(stdout);
        }
    }
    mpi_error=MPI_Finalize();
    if (MPI_SUCCESS!=mpi_error)
        return mpi_error;
    else
        return 0;
}